{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "335dc326",
   "metadata": {},
   "source": [
    "## 1.1.1 位置与分散程度的度量"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "bdda8cfc",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy.stats as st\n",
    "import pandas as pd\n",
    "import statsmodels.api as sm\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "plt.rcParams['font.sans-serif'] = ['SimHei']  # 用来正常显示中文标签\n",
    "plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "09c61415",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "62.36\n",
      "62.36\n",
      "学生体重的[10%, 20%, 40%, 60%, 80%, 100%]分位数： [52.76 56.98 62.2  64.   67.32 75.  ]\n",
      "体重数据方差的估计为：52.71，无偏估计为：56.47\n",
      "体重数据标准差的估计为：7.26，无偏估计为：7.51\n",
      "体重数据的变异系数为： 12.05 %\n",
      "体重数据的极差：27.60\n",
      "体重数据的标准误：1.94\n"
     ]
    }
   ],
   "source": [
    "# 示例：某学校15个学生体重（单位：公斤）抽样调查数据\n",
    "weights = np.array([75.0, 64.0, 47.4, 66.9, 62.2, 62.2, 58.7,\n",
    "                   63.5, 66.6, 64.0, 57.0, 69.0, 56.9, 50.0, 72.0])\n",
    "\n",
    "# 均值\n",
    "w_mean = np.mean(weights)\n",
    "w_mean2 = weights.mean()\n",
    "print(w_mean)\n",
    "print(w_mean2)\n",
    "\n",
    "# 限定范围内的数据求均值（截断 60 - 70）\n",
    "limitedMean = st.tmean(weights, (60, 70))\n",
    "\n",
    "sorted_weig = sorted(weights, reverse=True)  # reverse 的缺省值为False\n",
    "\n",
    "# 中位数，重要的统计值\n",
    "# 对称分布，比如T分布和正态分布，均值和中位数很接近，偏态分布的二者相差比较大，比如F分布\n",
    "median_weig = np.median(weights)\n",
    "\n",
    "# 分位数\n",
    "quantiles = np.quantile(weights, [0.1, 0.2, 0.4, 0.6, 0.8, 1])\n",
    "print('学生体重的[10%, 20%, 40%, 60%, 80%, 100%]分位数：', quantiles)\n",
    "\n",
    "# 方差、标准差、极差、标准误\n",
    "# 注意：方差与方差的无偏估计之间的计算区别\n",
    "v = np.var(weights) # 有偏估计或样本方差\n",
    "v_unb = st.tvar(weights) # 无偏估计\n",
    "print('体重数据方差的估计为：%0.2f，无偏估计为：%0.2f' % (v, v_unb))\n",
    "\n",
    "# 注意标准差与标准差的无偏估计之间的计算区别\n",
    "s = np.std(weights)  # 有偏估计或样本标准差\n",
    "s_unb = st.tstd(weights) # 无偏估计\n",
    "print('体重数据标准差的估计为：%0.2f，无偏估计为：%0.2f' % (s, s_unb))\n",
    "\n",
    "cv = s_unb / w_mean * 100  # 变异系数，无量纲，用百分数表示\n",
    "print('体重数据的变异系数为：', np.round(cv, 2), '%')\n",
    "\n",
    "# 极差与标准误\n",
    "R_weights = np.max(weights) - np.min(weights)  # 极差 = 最大值 - 最小值\n",
    "print('体重数据的极差：%0.2f' % R_weights)\n",
    "\n",
    "sm_weights = st.tstd(weights) / np.sqrt(len(weights))  # 标准误：数据标准差（无偏）/ 数据量 ** 0.5\n",
    "print('体重数据的标准误：%0.2f' % sm_weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "068d4c17",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "学生的身高、体重、胸围、坐高（前5个）：\n",
      " [[148  41  72  78]\n",
      " [139  34  71  76]\n",
      " [160  49  77  86]\n",
      " [149  36  67  79]\n",
      " [159  45  80  86]]\n",
      "\n",
      "学生的平均身高、平均体重、平均胸围、平均坐高分别为：\n",
      " 149.0, 38.7, 72.2, 79.4,\n"
     ]
    }
   ],
   "source": [
    "# 更复杂的例子：多维数据或矩阵求上述各统计量\n",
    "# 数据是从某学校抽样30个学生的身高、体重、胸围和坐高等抽样数据\n",
    "\n",
    "# 身高\n",
    "x1 = np.array([148, 139, 160, 149, 159, 142, 153, 150, 151, 139,\n",
    "              140, 161, 158, 140, 137, 152, 149, 145, 160, 156,\n",
    "              151, 147, 157, 147, 157, 151, 144, 141, 139, 148])\n",
    "\n",
    "# 体重\n",
    "x2 = np.array([41, 34, 49, 36, 45, 31, 43, 43, 42, 31,\n",
    "              29, 47, 49, 33, 31, 35, 47, 35, 47, 44,\n",
    "              42, 38, 39, 30, 48, 36, 36, 30, 32, 38])\n",
    "\n",
    "# 胸围\n",
    "x3 = np.array([72, 71, 77, 67, 80, 66, 76, 77, 77, 68,\n",
    "              64, 78, 78, 67, 66, 73, 82, 70, 74, 78,\n",
    "              73, 73, 68, 65, 80, 74, 68, 67, 68, 70])\n",
    "\n",
    "\n",
    "# 坐高\n",
    "x4 = np.array([78, 76, 86, 79, 86, 76, 83, 79, 80, 74,\n",
    "              74, 84, 83, 77, 73, 79, 79, 77, 87, 85,\n",
    "              82, 78, 80, 75, 88, 80, 76, 76, 73, 78])\n",
    "\n",
    "# 数据分析、统计建模、机器学习以及深度学习等应用领域一般将数据存储为列向量\n",
    "# numpy、scipy、pandas、statsmodels、sklearn、tensorflow、pytorch等基本上都是处理列向量\n",
    "# 当然也有例外，比如numpy求随机向量之间的协方差矩阵时，则是按照行向量进行计算的\n",
    "\n",
    "# 将x1,x2,x3,x4四个向量合并存储为矩阵，并转置为列向量\n",
    "stu_data = np.matrix([x1, x2, x3, x4]).T\n",
    "print('学生的身高、体重、胸围、坐高（前5个）：\\n', stu_data[0:5])\n",
    "\n",
    "# 注意：stu_data.mean(0)的用法，通过numpy的函数生成数据对象，numpy很多函数就会注入数据中\n",
    "# 此处直接调用数据对象的函数求平均值。函数的参数0表示列向量方向，如果为1则是行向量方向。\n",
    "# 也可以这样调用：np.mean(stu_data, axis=0)，结果是一样的\n",
    "# 其他统计量计算类似\n",
    "\n",
    "# ravel 对数据展平\n",
    "stu_mean = np.round(stu_data.mean(0), 1).ravel() # 将二维矩阵展平为一维向量，数据类型转变为numpy数据\n",
    "print('\\n学生的平均身高、平均体重、平均胸围、平均坐高分别为：\\n %.1f, %.1f, %.1f, %.1f,' % (stu_mean[0], stu_mean[1], stu_mean[2], stu_mean[3]))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a10cfe6",
   "metadata": {},
   "source": [
    "## 1.1.2 关系度量"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "4907a6bb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "学生身高、体重、胸围和坐高之间的协方差与相关系数矩阵分别如下：\n",
      "\n",
      " [[53.52 40.79 27.59 28.76]\n",
      " [40.79 41.73 29.83 24.36]\n",
      " [27.59 29.83 26.53 17.22]\n",
      " [28.76 24.36 17.22 18.24]] \n",
      "\n",
      " [[1.   0.86 0.73 0.92]\n",
      " [0.86 1.   0.9  0.88]\n",
      " [0.73 0.9  1.   0.78]\n",
      " [0.92 0.88 0.78 1.  ]]\n"
     ]
    }
   ],
   "source": [
    "## 协方差矩阵：covariance\n",
    "cov_stu = np.cov(stu_data.T) # 方差-协方差矩阵\n",
    "\n",
    "## 相关系数矩阵：correlation coefficient\n",
    "rou_stu = np.corrcoef(stu_data.T)\n",
    "\n",
    "# 四个随机向量之间的相关程度很高，尤其是身高和坐高之间的相关性最高\n",
    "print('学生身高、体重、胸围和坐高之间的协方差与相关系数矩阵分别如下：\\n\\n', np.round(cov_stu, 2), '\\n\\n', np.round(rou_stu, 2))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf3e257f",
   "metadata": {},
   "source": [
    "## 1.1.3 分布形态的度量\n",
    "\n",
    "#### 偏度计算：\n",
    "- 1. 偏度表示曲线是向左偏或右偏，又称为正偏态或负偏态\n",
    "- 2. 偏度越接近0，岳父和正态分布的曲线。\n",
    "- 3. 偏度小于0称分布具有负偏态，也称左偏态；反正就是正偏态或右偏态"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "c5d3d9b7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Pandas计算公式 手工计算以及调用函数计算结果：\n",
      "skew1: -0.4299560852122288 , skew_pandas: -0.429956085212229\n",
      "\n",
      "Scipy 计算公式手工计算以及调用函数计算结果（无修正）：\n",
      "skew2: -0.38570725236501335 , skew_scipy: -0.38570725236501346\n",
      "\n",
      " Scypy进行修正后的偏度： -0.429956085212229\n"
     ]
    }
   ],
   "source": [
    "# 偏度计算公式\n",
    "n = len(weights)\n",
    "\n",
    "# 三阶矩，其他各阶矩的计算以此类推\n",
    "u3 = np.sum((weights - w_mean) ** 3) / n\n",
    "# 四阶矩\n",
    "u4 = np.sum((weights - w_mean) ** 4) / n\n",
    "\n",
    "# 使用总体标准差的无偏估计，计算的偏度是修正后的偏度\n",
    "skew1 = ((n**2)*u3)/((n-1)*(n-2)*(s_unb**3))\n",
    "\n",
    "# pandas 计算是修正后偏度\n",
    "pd_weights = pd.Series(weights)\n",
    "skew_pandas = pd_weights.skew()\n",
    "print('Pandas计算公式 手工计算以及调用函数计算结果：')\n",
    "print('skew1:', skew1, ', skew_pandas:', skew_pandas)\n",
    "\n",
    "# 无修正偏度的手工计算，使用样本标准差\n",
    "skew2 = np.sum((weights - w_mean) ** 3) / ((s ** 3) * n)\n",
    "\n",
    "# scipy计算公式和结果\n",
    "print('\\nScipy 计算公式手工计算以及调用函数计算结果（无修正）：')\n",
    "skew_scipy = st.skew(weights)\n",
    "print('skew2:', skew2, ', skew_scipy:', skew_scipy)\n",
    "\n",
    "# 1. 使用Scipy的skew函数，如果将第二个参数bias设置为False，计算结果就和Pandas完全相同了。\n",
    "#    bias参数表示是否修正，如果为False表示修正，反正则不修正。\n",
    "# 2. 总体上感觉修正后偏度比较准确，但是很多场合仍用无修正的偏度进行统计量的计算。\n",
    "# 3. StatsModels的线性回归模型对残差的正态分布性（Jarque-Bera、Omnibus检验等）\n",
    "#    进行检验时，使用的偏度就是无修正的，包括峰度也是无修正的。\n",
    "\n",
    "skew_scipy_bias = st.skew(weights, bias=False)\n",
    "print('\\n Scypy进行修正后的偏度：', skew_scipy_bias)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "427b0b04",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      " Pandas计算峰度： 0.09653947135209195 \n",
      "\n",
      " Scipy计算峰度（修正后：） 0.09653947135209329 \n",
      "\n",
      " Scipy计算峰度（无修正：） -0.3077671538797926\n"
     ]
    }
   ],
   "source": [
    "# 峰度的计算：\n",
    "# 1. 峰度表示曲线是扁平态（低峰态）还是尖峰态\n",
    "# 2. 正常值有两种定义：Fisher定义该值为0；Pearson定义为3\n",
    "# 3. 按照Fisher定义，峰度=0表示正好符合正态分布的曲线；大于0表示峰比较尖，反之表示比较平。\n",
    "\n",
    "\n",
    "# 峰度计算，StatsModels多使用无修正的峰度\n",
    "# 手工实现留作练习\n",
    "kurt_pandas = pd_weights.kurt()\n",
    "kurt_scipy = st.kurtosis(weights, bias=False)\n",
    "kurt_scipy_bias = st.kurtosis(weights, bias=True) # True是bias的缺省值\n",
    "print('\\n Pandas计算峰度：', kurt_pandas, '\\n\\n Scipy计算峰度（修正后：）', kurt_scipy, '\\n\\n Scipy计算峰度（无修正：）', kurt_scipy_bias)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e790bc45",
   "metadata": {},
   "source": [
    "## 1.1.4 数据特性的总括\n",
    "\n",
    "#### 数据总括：\n",
    "- 1. 数据特性总括，包括：最小最大值、均值、方差、偏度、峰度。\n",
    "- 2. 此处使用了修正选项，偏度和峰度都是修正后的值\n",
    "- 3. 数据分布的正态性检验与分布拟合检验"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "520eb6ac",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "学生体重数据的总括描述： DescribeResult(nobs=15, minmax=(47.4, 75.0), mean=62.36, variance=56.472571428571435, skewness=-0.429956085212229, kurtosis=0.09653947135209329)\n",
      "\n",
      " 学生体重数据的正态性检验： ShapiroResult(statistic=0.9686222076416016, pvalue=0.8371511101722717)\n",
      "\n",
      " 测试不服从正态分布的数据： ShapiroResult(statistic=0.5549201965332031, pvalue=0.0001440312626073137)\n"
     ]
    }
   ],
   "source": [
    "print('学生体重数据的总括描述：', st.describe(weights, bias=False))\n",
    "\n",
    "# 正态性检验：主要考察p值，如果p值>0.05，则不能拒绝原假设，即数据服从正态分布\n",
    "# 此处的p值=0.837，远大于0.05，说明体重数据服从正态分布性\n",
    "# 最常用的函数是夏皮罗（shapiro）函数检验\n",
    "# Scipy包含多种正态分布检验的函数，后面课程会陆续接触到。\n",
    "\n",
    "print('\\n 学生体重数据的正态性检验：', st.shapiro(weights))\n",
    "print('\\n 测试不服从正态分布的数据：', st.shapiro([1,2,3,4,900]))  # p值远小于0.05，拒绝原假设"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "179e9c96",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "检验数据是否服从某种分布： KstestResult(statistic=0.5019030078083166, pvalue=3.532570881297376e-12)\n",
      "\n",
      "检验数据是否服从某种分布： KstestResult(statistic=0.09418170098517131, pvalue=0.7313575462160135)\n"
     ]
    }
   ],
   "source": [
    "# 经验分布的检验方法：Kolmogorov-Smirnov检验法\n",
    "# 该函数检验数据是否服从某个类型的概率分布函数，不仅仅是正态分布\n",
    "# 在数据分析和机器学习领域，很多算法的前提要求数据服从正态分布\n",
    "# kstest函数原假设两个独立样本数据来自同一个连续分布\n",
    "# 生成一个服从F分布的测试数据：\n",
    "# 1. 第一个函数测试其是否服从自由度为3的t分布；\n",
    "# 2. 第二个函数测试是否服从自由度为2，9的F分布。\n",
    "\n",
    "# 生成服从F分布，自由度为（2，9）的随机数据\n",
    "f_data = st.f.rvs(size=50, dfn=2, dfd=9)\n",
    "\n",
    "# 检验上述数据是否服从自由度为3的T分布，结果很显然拒绝服从该分布的原假设\n",
    "print('\\n检验数据是否服从某种分布：', st.kstest(f_data, 't', (3,))) # (3,)表示T分布的自由度\n",
    "print('\\n检验数据是否服从某种分布：', st.kstest(f_data, 'f', (2, 9))) # 接受原假设服从自由度为2，9的F分布"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "09cd740b",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
